Transcriptome analyses describe the consequences of persistent HIF-1 over-activation in Caenorhabditis elegans

Metazoan animals rely on oxygen for survival, but during normal development and homeostasis, animals are often challenged by hypoxia (low oxygen). In metazoans, many of the critical hypoxia responses are mediated by the evolutionarily conserved hypoxia-inducible transcription factors (HIFs). The stability and activity of HIF complexes are strictly regulated. In the model organism C. elegans, HIF-1 stability and activity are negatively regulated by VHL-1, EGL-9, RHY-1 and SWAN-1. Importantly, C. elegans mutants carrying strong loss-of-function mutations in these genes are viable, and this provides opportunities to interrogate the molecular consequences of persistent HIF-1 over-activation. We find that the genome-wide gene expression patterns are compellingly similar in these mutants, supporting models in which RHY-1, VHL-1 and EGL-9 function in common pathway(s) to regulate HIF-1 activity. These studies illuminate the diversified biological roles played by HIF-1, including metabolism and stress response. Genes regulated by persistent HIF-1 over-activation overlap with genes responsive to pathogens, and they overlap with genes regulated by DAF-16. As crucial stress regulators, HIF-1 and DAF-16 converge on key stress-responsive genes and function synergistically to enable hypoxia survival.


Introduction
As the electron acceptor during oxidative phosphorylation for energy production, oxygen is vital to all the aerobic organisms.Insufficient oxygen availability not only decreases energy production but also changes the redox environment for cellular biochemical reactions [1,2].During normal development and disease states, organisms are often challenged by hypoxia.In metazoans, HIFs have evolutionarily conserved roles in mediating critical transcriptional responses to hypoxia.In mammals, HIF complexes regulate genes functioning in angiogenesis, erythropoiesis and glycolysis to assist hypoxia survival [3][4][5][6].HIF's function and regulation have potential therapeutic significance in hypoxia-related diseases, such as cancer and stroke.HIF complexes are heterodimers composed of an α subunit and a β subunit, and both subunits are bHLH (basic-helix-loop-helix)-PAS (PER/ARNT/SIM) domain proteins [7].The human genome encodes three HIFα subunits: HIF-1, -2 and -3α, and three HIFβ/ARNT subunits: HIF-1, -2 and -3β [8].While HIFβ has multiple bHLH-PAS dimerization partners and is relatively stable and abundant, HIFα is short-lived under well-oxygenated conditions.Thus, HIFα is dedicated to regulate the expression of oxygen-sensitive genes [9][10][11].The oxygen-dependent HIFα degradation pathway is conserved and well established.When oxygen levels are high, prolyl hydroxylase domain proteins (PHDs) hydroxylate specific proline residues on HIFα, using oxygen as co-substrate.The hydroxylated HIFα is targeted for proteasomal degradation by an E3 ubiquitin ligase containing the tumor suppressor von Hippel-Lindau (VHL) [2,10].
The PHD-VHL system for HIF stability regulation is conserved and simplified in C. elegans, too.While humans have three PHDs, C. elegans has only one counterpart encoded by egl-9; and the single VHL homolog in C. elegans is encoded by vhl-1 [14].
These mutants provide an opportunity to employ multiple genetic backgrounds to overactivate HIF-1 and determine the downstream effects.This also provides some insights to the molecular networks that enable animals to respond to diverse stresses.We answer these questions by comparing the genome-wide transcriptional profiles in these HIF-1 negative regulator mutants.
the overlaps are statistically significant (p-values <2.2E-16 by Fisher's exact tests) (Fig 1A).In sum, the genes that were up-regulated in these three mutants overlapped extraordinarily to illuminate the consequences of long-term HIF-1 over-activation.
In sum, the microarray data revealed that the gene expression patterns in the three HIF-1 high-activity mutants were strikingly similar.This observation supported existing models that RHY-1, VHL-1 and EGL-9 functioned in common pathway(s) to regulate HIF-1 activity.

Genes regulated by HIF-1 in the HIF-1 negative regulator mutants and under hypoxia
We anticipated some overlap between genes that were mis-regulated in egl-9, vhl-1 or rhy-1 loss-of-function mutants and genes that had been shown to be induced by hypoxia in a hif-1dependent manner.To test this hypothesis, we asked whether the 104 genes that were up-regulated in all 3 mutants (S11 Table) included genes that had been shown to be positively regulated by HIF-1 in hypoxic conditions [85].The overlap is striking, as illustrated in S1 Fig, and includes genes for lipid metabolism (ZK550.6 and gbh-2), for propionic acid metabolism (mce-1 and mmcm-1), H 2 S and HCN detoxification (cysl-2, ethe-1 and sqrd-1), gluconeogenesis (pck-1), and protein synthesis regulation (efk-1), also for collagen synthesis (phy-2) (S13 Table ).We also looked for overlaps between the 95 genes that were expressed at lower levels in all three mutants (S13 Table) and genes that had been shown to be repressed by hypoxia in a hif-1-dependnet manner [85].There were fewer genes that were common to all of these gene sets (S1 Fig), and they included T28A11.2 (hypothetical protein), acdh-2 (Acyl CoA dehydrogenase), and acs-2 (fatty acid CoA synthetase family) (S14 Table ).

Consequences of persistent HIF-1 over-activation
To more fully understand the consequences of persistent HIF-1 over-activation, we explored the functions the 104 genes that were commonly up-regulated in vhl-1(ok161), egl-9(sa307) and rhy-1(ok1402) (S11 Table ), and the 95 genes that were commonly down-regulated in these three HIF-1 high-activity mutants (S12 Table ).The enriched biological terms for genes commonly up-regulated in these three HIF-1 high-activity mutants identified by WormCat at Category 1 level were stress response and metabolism (Table 1).The enriched biological term for these genes at Category 2 level was stress response: pathogen (Table 1).And the enriched biological term for them at Category 3 level was stress response: pathogen: CUB.Genes enriched in these functional categories are listed in Table 1.For genes commonly down-regulated in the three HIF-1 high activity mutants, there was no enriched biological term at Category 1 level with Bonferroni false discovery rate setting at 0.01.The enriched biological term for these genes at Category 2 level were: signaling: hedgehog-like and metabolism: lipid (Table 2).And the enriched biological term for them at Category 3 level were: signaling: hedgehog-like and unassigned: regulated by multiple stresses (Table 2).Genes enriched in each functional category are listed in Table 2.

Testing whether nhr-49 hypoxia pathway functions downstream of vhl-1, egl-9 or rhy-1
A recent study has demonstrated that the transcription factor NHR-49 acts in parallel to HIF-1 to regulate hypoxia response, and 83 genes were identified as NHR-49-dependent and HIF-1-independent hypoxia inducible genes [89].We asked whether these 83 NHR-49 targets were differentially expressed in vhl-1(ok161), rhy-1(ok1402) and egl-9(sa307) mutants.We found that none of these overlaps were significant (Table 3).This is consistent with models in which HIF-1 and its upstream regulators act in parallel to NHR-49.
We compared genes regulated by DAF-16 [101] with those regulated by HIF-1 over-activation in vhl-1(ok161), swan-1(ok267);vhl-1(ok161), egl-9(sa307) or rhy-1(ok1402) to see whether they had significant overlaps.There was significant overlap between genes up-regulated in each of these four HIF-1 negative regulator mutants and genes that have been shown to be positively regulated by DAF-16 (Table 5).But there were no significant overlaps between genes down-regulated in these four mutants and genes down-regulated by DAF-16 (Table 5).
We further identified 18 genes that were commonly regulated by DAF-16 and in all the four HIF-1 high-activity mutants.Among them, 15 genes were up-regulated by DAF-16 and up-regulated in all the four HIF-1 high-activity mutants; this overlap is significant (pvalue = 7.26E-12, by Fisher's exact test).Three genes were down-regulated by DAF-16 and down-regulated in all the four HIF-1 high-activity mutants; this overlap is not significant (pvalue = 0.11, by Fisher's exact test).The expression profiles of these 18 genes in vhl-1(ok161), swan-1(ok267);vhl-1(ok161), egl-9(sa307) and rhy-1(ok1402) are illustrated in a heat map in We further investigated the requirement for daf-16 in moderate hypoxia adaptation (0.5% oxygen, 21˚C).Consistent with prior studies, we found that 100% of N2 eggs hatched within 24 hours in hypoxia, and 100% developed to adulthood within 72 hours.By contrast, only 52%  of hif-1(ia04) mutant eggs hatched within 24 hours in hypoxia, and only 10% developed to adulthood within 72 hours.As shown in Table 6, loss-of-function mutations in daf-16 also impaired embryogenesis and larval development in hypoxia.For the two daf-16 loss-of-function alleles (daf-16(mu86) and daf-16 (mgDf50)) tested, 43-46% of daf-16-deficient eggs hatched within 24 hours in hypoxia, and 30-39% developed to adulthood within 72 hours (Table 6).The double mutants of daf-16 and hif-1 (daf-16(mgDf50);hif-1(ia04)) were most severely affected by hypoxia, only 8% of eggs hatched within 24 hours in hypoxia, and only 0.4% developed to adulthood within 72 hours (Table 6).Note, in room air, the hatched or adult rates for all of these mutant genotypes were the same as those in N2 wild type: 100% of the eggs hatched within 24 hours, and 100% of them developed to adulthood within 72 hours (Table 6).Using the standard provided by Kirienko et al. [104], these data suggest that HIF-1 and DAF-16 acted synergistically to protect C. elegans from hypoxia insult.

Discussion
Prior studies have shown that RHY-1, EGL-9 and VHL-1 regulate HIF-1 activity [14,39,50,83], but the extent to which their functions overlap was not fully understood.Here, by employing whole-genome transcriptome analyses, we are able to shed more light on this important regulatory network.The analyses of genes that were differentially expressed in rhy-1, egl-9 and vhl-1 mutants also extend our understanding of the consequences of multi-generational HIF-1 activation and provides insights to how HIF-1 pathway interacts with other pathways that mediate stress responses.

The pleiotropic consequences of persistent HIF-1 over-activation account for its diversified biological roles
Misregulation of HIF-1 disrupts multiple facets of C. elegans' life, including hypoxia adaptation, hydrogen sulfide and hydrogen cyanide resistance, heat resistance, heavy metal toxicity tolerance, pathogen response, egg-laying defect, brood size and aging [22, 29, 30, 34-36, 38, 39, 44, 50, 61-64, 81, 82].Our analyses provide a more complete understanding of gene expression changes that underpin these diverse phenotypes.The common set of genes misregulated in vhl-1(ok161), egl-9(sa307) and rhy-1(ok1402) are involved in the metabolism of sugars, lipids, amino acids, and sulfur.Consistent with this, it have been shown that the metabolomic profiles of various amino acids, carbohydrates, lipids and nucleotides are changed in egl-9(sa307) mutants [106].They also have roles in hypoxia and innate immune response.These various biological processes regulated by persistent HIF-1 over-activation suggest models for why HIF-1 can play a variety of biological roles besides regulating hypoxia response.For example, the up-regulated sulfur metabolism by HIF-1 over-activation is relevant to HIF-1's roles in hydrogen sulfide and hydrogen cyanide resistance [29,30,[32][33][34].Also, the up-regulated innate immune response by HIF-1 over-activation explains its roles in pathogen resistance [35, 36, 38-41, 81, 104], which is further supported by the overlaps of genes regulated by HIF-1 over-activation with genes responsive to toxic Cry5B, P. aeruginosa PA14, S. aureus and Y. pestis [36,92,94,95], and with genes regulated by PMK-1 and SEK-1, the major players for pathogenic P. aeruginosa PA14 resistance [91].We recognized that the methods we used to identify the overlaps between differentially expressed gene lists has its own limitation.When doing genome-wide analyses, researchers rely on statistical cut offs to focus inquiry on the genes for which the statistical data is the strongest.Additional repetitions would no doubt identify more genes as differentially expressed.This dataset provides a starting point for further study.

HIF-1 interacts with DAF-16 to promote hypoxia adaptation
DAF-16 and HIF-1 are both important stress regulators in C. elegans.In this study, we demonstrate that DAF-16 and HIF-1 have complicated complementary and overlapping roles in stress response.On one hand, they converge on key stress-responsive genes.On the other hand, they function synergistically to promote hypoxia adaptation.This suggests that DAF-16 is largely required outside the HIF-1 pathway to promote hypoxia resistance.Intriguingly, the mRNA levels of daf-18, the positive regulator of DAF-16, are up-regulated in swan-1(ok267); vhl-1(ok161), egl-9(sa307) and rhy-1(ok1402) mutants.This suggests a mechanism for the crosstalk between HIF-1 and DAF-16: HIF-1 may enhance DAF-16's function by increasing the expression of daf-18.While beyond the scope of the experiments reported here, future studies might explore the ways in which these two important stress-response pathways regulate each other and converge to enable the organism to adapt and survive adverse conditions.

Gene expression microarray experiment
Randomized complete block design was followed for the microarray experiment, with three biological replicates treated as three blocks.Each block included eight treatments: N2 wild type, N2 wild type with hypoxia treatment, hif-1(ia04) loss-of-function mutants, hif-1(ia04) loss-of-function mutants with hypoxia treatment, vhl-1(ok161) loss-of-function mutants, rhy-1 (ok1402) loss-of-function mutants, egl-9(sa307) loss-of-function mutants and swan-1(ok267); vhl-1(ok161) loss-of-function double mutants.For each treatment, about 1,000 synchronized L4-stage larvae were pooled as one experimental unit to get sufficient RNA for hybridization.Total RNA isolation was performed using Trizol (Invitrogen) and RNeasy Mini Kit (Qiagen).RNA quality was checked with an Agilent 2100 BioAnalyzer (Agilent Technologies).The RNA integrity numbers (RINs) for all the samples used in this study were greater than 9.0.The total RNA isolated from one experimental unit was hybridized onto one Affymetrix GeneChip1 C. elegans Genome array (Affymetrix, part number 900383).Probe synthesis, labeling, hybridization, washing, staining and scanning were performed by the GeneChip facility at Iowa State University.In brief, the total RNA was synthesized to biotin-labeled aRNA using the Gene-Chip1 3' IVT Express Kit (Affymetrix, part number 901229) and hybridized to the array.The arrays were washed and stained in the GeneChip1 fludics station 450 and scanned with Gene-Chip1 scanner 3000 7G.The Affymetrix1 GeneChip1 Command Console™ (AGCC) software was used to generate probe cell intensity data (.CEL) files.The resulting CEL files were normalized and summarized using the robust multichip average (RMA) algorithm [108] in R package (R Core Team, Vienna, Austria, 2016).An analysis of variance (ANOVA) model was then fitted to the summarized expression measures, with the block (three levels) and the treatment (eight levels) treated as fixed effect factors following the experimental design.Residual model diagnostics identified no severe violations of the model assumptions.Linear contrasts of treatment means were tested using the general F-test.To account for multiplicities of hypothesis testing, conservative estimates of false discovery rates (FDRs) were calculated according to the q-value procedure of Storey and Tibshirani [109].Differentially expressed probesets were defined as q-value � 0.05 and fold change � 1.6.Probesets were converted to genes using the Affymetrix annotation file "Celegans.na36.annot.csv".To deal with redundancy and count the number of unique genes detected on the array, we kept one probeset per gene and one gene per probeset.In this way, the total number of unique genes detected on the array was 18, 011.For the purpose of reference, the original complete lists of gene(s) annotated to each probeset were kept in S1-S8 Tables.The complete analysis results for all the conditions (N2 wild type, N2 wild type with hypoxia treatment, hif-1(ia04) loss-of-function mutants, hif-1 (ia04) loss-of-function mutants with hypoxia treatment, vhl-1(ok161) loss-of-function mutants, rhy-1(ok1402) loss-of-function mutants, egl-9(sa307) loss-of-function mutants and swan-1(ok267);vhl-1(ok161) loss-of-function double mutants) and all the probesets on the microarray are provided in S1 Table .Gene expression changes in N2 wild type animals and hif-1(ia04) mutants under hypoxia and room air, and the comparison of gene expressions under hypoxia and in the HIF-1 negative regulator mutants, also the expression of HIF-1 direct targets identified by CHIP-seq have been presented in a related study [85].This study focuses on the genes expression changes in vhl-1(ok161), egl-9(sa307), rhy-1(ok1402) and swan-1(ok267);vhl-1(ok161) double mutants.The microarray raw and probeset summary data had been deposited to NCBI's Gene Expression Omnibus, the accession number was GSE228851.

Gene function annotation and enrichment analyses
WormCat online tool (www.wormcat.com)[110] was used to annotate the enriched biological terms associated with microarray-selected genes.The enriched biological terms were at Bonferroni false discovery rate cut off of 0.01.

Heat maps
Heat maps for gene expression profiles were generated using the PermutMatrix graphical analysis program [111,112].Average linkage clustering was performed with the fold changes compared to N2. Green color represented negative values, and red color represented positive values.The color intensities corresponded to the magnitudes of fold changes.Other parameters were set as default.

Gene lists overlap testing
Fisher's exact test was performed to test whether the overlap between two gene lists was significant or not.The total number of 18, 011 genes detected on the microarray was used as the population size.The significant overlap is at p-value < 0.001.

Hypoxia development and survival assays
For each genotype, the room air and hypoxia treatments were performed in parallel at 21˚C.For each treatment, 20 young adults (one day after L4 molt) were used as parents to lay eggs on one NGM plate seeded with OP50 for 30 minutes.After counting the eggs laid, the plates were kept in room air or put into a sealed plexiglass chamber with constant hypoxic gas flow.Compressed air and 100% nitrogen were mixed to achieve 0.5% oxygen, and gas flow was controlled by an oxygen sensor [84].After 24 hours, the un-hatched eggs were counted for both treatments.After that, the plates for both treatments were maintained in room air.The adult worms were counted 72 hours after the eggs had been laid.The data collection time points were set to match the development rate of N2 eggs in room air: they hatch within 24 hours and reach adulthood within 72 hours.The experiments were performed with three biological replicates.To test the effect of hypoxia on animal development and survival, the binary hatched vs. un-hatched or adult vs. non-adult data were analyzed for each genotype by fitting a generalized linear model using a logit link function with JMP 9 statistical software (SAS Institute Inc., Cary, NC, 2010).The replicate (three levels) and the treatment (two levels) were used as factors in the model.

Table 6 . hif-1 and daf-16 functioned synergistically to protect C. elegans in hypoxia.
n is the total number of animals assayed from three biological replicates.